/* fc_layer.c */
#include "fc_layer.h"
#include "arm_math.h" /* cmsis-dsp. */
#include <math.h> /* double exp(double); */
#include <stdlib.h> /* rand(x) */

#define FC_LAYER_VECTOR_LEN_MAX   256u

/*
 * fc_layer_tmp1_arr内存是给最底层的基础计算函数使用的, 包括: arm_tanh_f32()和arm_dtanh_f32()
 * 调用基础计算函数的高层函数不能使用, 否则其暂存数据会被基础函数的执行而被充掉.
 */
float fc_layer_tmp1_arr[FC_LAYER_VECTOR_LEN_MAX];
float fc_layer_tmp2_arr[FC_LAYER_VECTOR_LEN_MAX];
float fc_layer_tmp3_arr[FC_LAYER_VECTOR_LEN_MAX];

/* exp(x) = e^x. */
void arm_exp_f32(float *input, float *output, uint32_t len)
{
    for (uint32_t i = 0u; i < len; i++)
    {
        *output = exp(*input);
        output++;
        input++;
    }
}

/* div(x1, x2) = x1 / x2. */
void arm_div_f32(float *input1, float *input2, float *output, uint32_t len)
{
    for (uint32_t i = 0u; i < len; i++)
    {
        *output = (*input1) / (*input2);
        output++;
        input1++;
        input2++;
    }
}

/* tanh(x) = (exp(x)-exp(-x)) / (exp(x)+exp(-x)). */
void arm_tanh_f32(float *input, float *output, uint32_t len)
{
    /* tmp1 = exp(x). */
    arm_exp_f32(input, fc_layer_tmp1_arr, len);

    /* tmp2 = exp(-x). */
    arm_negate_f32(input, fc_layer_tmp2_arr, len);
    arm_exp_f32(fc_layer_tmp2_arr, fc_layer_tmp2_arr, len);

    /* output = tmp1 - tmp2. */
    arm_sub_f32(fc_layer_tmp1_arr, fc_layer_tmp2_arr, output, len);

    /* tmp3 = tmp1 + tmp2. */
    arm_add_f32(fc_layer_tmp1_arr, fc_layer_tmp2_arr, fc_layer_tmp3_arr, len);

    /* output = output / tmp4. */
    arm_div_f32(output, fc_layer_tmp3_arr, output, len);
}

/* dtanh(x) = 1.0f - x^2. */
void arm_dtanh_f32(float *input, float *output, uint32_t len)
{
    /* tmp1 = input * input. */
    arm_mult_f32(input, input, fc_layer_tmp1_arr, len);

    /* output = 1.0. */
    arm_fill_f32(1.0f, output, len);

    /* output = output - tmp1. */
    arm_sub_f32(output, fc_layer_tmp1_arr, output, len);
}

void arm_rand_f32(float *output, float min, float max, uint32_t len)
{
    //return ((float(rand()) / float(RAND_MAX)) * (Max - Min)) + Min;
    for (uint32_t i = 0u; i < len; i++)
    {
        *output = ((float)rand() / (float)RAND_MAX) * (max - min) + min;
        output++;
    }
}

/******************************************************************************
 * fc_layer的功能函数.
 *****************************************************************************/
void fc_layer_new_infer(fc_layer_t *layer, uint32_t dim_in, uint32_t dim_out, float *W, float *B)
{
    layer->dim_in = dim_in;
    layer->dim_out = dim_out;
    layer->W = W;
    layer->B = B;
}

void fc_layer_infer(fc_layer_t *layer, float *in_val, float *val_out)
{
    /* 1. 输入乘以权重得到初级输出. */
    float *pW = layer->W;
    //float *pO = val_out;
    for (uint32_t i = 0u; i < layer->dim_out; i++)
    {
        arm_dot_prod_f32(in_val, pW, layer->dim_in, &val_out[i]); /* 先乘以权重. */
        //pO++;
        pW += layer->dim_in;
    }
    //pO = val_out;

    /* 2. 加上偏置得到二级输出. */
    arm_add_f32(val_out, layer->B, val_out, layer->dim_out);

    /* 3. 传入激活函数得到最终输出. */
    arm_tanh_f32(val_out, val_out, layer->dim_out);
}

bool fc_layer_new_backprop(fc_layer_t *layer, uint32_t dim_in, uint32_t dim_out, float *dW, float *dB)
{
    if ( (layer->dim_in != dim_in) || (layer->dim_out != dim_out) )
    {
        return false;
    }

    layer->dW = dW;
    layer->dB = dB;
    return true;
}

void fc_layer_clear_gradients(fc_layer_t *layer)
{
    arm_fill_f32(0.0f, layer->dW, layer->dim_out * layer->dim_in);
    arm_fill_f32(0.0f, layer->dB, layer->dim_out);
}

/*
 * 反向传播误差函数.
 * prv_out: 输入, 本层在之前推导得到的输出
 * prv_err: 输入, 本层预测的误差, 后一级向本级的反馈
 * in_val : 输入, 本级输入, 前级的输出
 * in_err : 输出, 向前级输出反馈的偏差(通过激活函数的梯度反向传播)

 * 注意, 本函数看到的prv_err不是用prv_out和真实值直接差出来的, 在nn的中间层没有真实值, 只有算到最后才能有真实值.
 * "真实的数据不是根据你网络设计出来的, 人家有自己的规律, 只是告诉你开始和结束. 你要自己找规律"
 * 此处的pre_err是从后一级反向传播反馈回来的, 整个网络从最后的输出级向前级逐层传播, 一边修正.
 * 本函数处理的数据涉及后级, 本级, 前级. 连接相邻的三层, 从后向前传播.
 */
float fc_layer_grad_in_arr[FC_LAYER_VECTOR_LEN_MAX]; /* 暂存反向传播导数. */
void fc_layer_backprop(fc_layer_t *layer, float *prv_out, float *prv_err, float *in_val, float *in_err)
{
    if (in_err != NULL)
    {
        arm_fill_f32(0.0f, in_err, layer->dim_in);
    }

    /* 得到反向传播导数的值(系数). dO * dtanh(y) */
    arm_dtanh_f32(prv_out, fc_layer_grad_in_arr, layer->dim_out);
    arm_mult_f32(prv_err, fc_layer_grad_in_arr, fc_layer_grad_in_arr, layer->dim_out);

    /* 此后, 将反向传播系数转化成权重和偏置的调整值. */

    /* 调整权重, 并且累计对输入的反馈. */
    for (uint32_t i_out = 0u; i_out < layer->dim_out; i_out++)
    {
        /* 每个神经元的输出将自己的偏差反馈到与自己相连的所有输入通路上. */
        for (uint32_t i_in = 0u; i_in < layer->dim_in; i_in++)
        {
            layer->dW[i_in * layer->dim_out + i_out] -= in_val[i_in] * fc_layer_grad_in_arr[i_out]; /* 对权重进行反馈. 以权重为变量求偏导数.*/
            if (in_err != NULL) /* 反向传播误差. */
            {
                in_err[i_in] -= layer->W[i_in * layer->dim_out + i_out] * fc_layer_grad_in_arr[i_out]; /* 对输入进行反馈. 以输入为变量求偏导数. */
            }
        }
    }

    /* 调整偏置. */
    arm_sub_f32(layer->dB, fc_layer_grad_in_arr, layer->dB, layer->dim_out);
}

void fc_layer_update(fc_layer_t *layer, float learn_rate)
{
    arm_scale_f32(layer->dW, learn_rate, layer->dW, layer->dim_out * layer->dim_in);
    arm_scale_f32(layer->dB, learn_rate, layer->dB, layer->dim_out);

    arm_sub_f32(layer->W, layer->dW, layer->W, layer->dim_out * layer->dim_in);
    arm_sub_f32(layer->B, layer->dB, layer->B, layer->dim_out);

    /* 一旦被更新到权重矩阵和偏置向量中去, 准备的矩阵调整量和偏置调整量就没用了.
     * 注意, 在更新之前, 矩阵调整量和偏置调整量是可以累积的. 可以积累多次之后再进行调整.
     */
    fc_layer_clear_gradients(layer);
}

/* EOF. */

